Jan-Philipp Kolb
27 April 2017
R Welt
Aktivität Nutzer
Das kann diese Schulung vermitteln:
Das kann sie nicht leisten:
Freie Beteiligung - modularer Aufbau (immer mehr Erweiterungspakete)
Der Download ist auf dieser Seite möglich:
Aber die meisten Menschen nutzen einen Editor oder ein graphical user interface (GUI).
Aus den folgenden Gründen:
Auf github sind alle Unterlagen für diesen Kurs zu finden.
Sechs Gründe Rstudio zu nutzen.
Wie man Rstudio nutzen kann.
date() und die R-version mit sessionInfo() heraus.date()## [1] "Thu Apr 27 22:45:09 2017"
sessionInfo()## R version 3.3.3 (2017-03-06)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 7 x64 (build 7601) Service Pack 1
##
## locale:
## [1] LC_COLLATE=German_Germany.1252 LC_CTYPE=German_Germany.1252
## [3] LC_MONETARY=German_Germany.1252 LC_NUMERIC=C
## [5] LC_TIME=German_Germany.1252
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## loaded via a namespace (and not attached):
## [1] backports_1.0.5 magrittr_1.5 rprojroot_1.2 tools_3.3.3
## [5] htmltools_0.3.5 yaml_2.1.13 Rcpp_0.12.6 stringi_1.1.1
## [9] rmarkdown_1.4 knitr_1.15.1 stringr_1.2.0 digest_0.6.12
## [13] evaluate_0.10
Vektoren und Zuweisungen
<- ist der Zuweisungsoperatorb <- c(1,2) # erzeugt ein Objekt mit den Zahlen 1 und 2mean(b) # berechnet den Mittelwert## [1] 1.5
Mit den folgenden Funktionen können wir etwas über die Eigenschaften des Objekts lernen:
length(b) # b hat die Länge 2## [1] 2
str(b) # b ist ein numerischer Vektor## num [1:2] 1 2
| Funktion | Bedeutung | Beispiel |
|---|---|---|
| length() | Länge | length(b) |
| max() | Maximum | max(b) |
| min() | Minimum | min(b) |
| sd() | Standardabweichung | sd(b) |
| var() | Varianz | var(b) |
| mean() | Mittelwert | mean(b) |
| median() | Median | median(b) |
Diese Funktionen brauchen nur ein Argument.
Andere Funktionen brauchen mehr:
| Argument | Bedeutung | Beispiel |
|---|---|---|
| quantile() | 90 % Quantile | quantile(b,.9) |
| sample() | Stichprobe ziehen | sample(b,1) |
max(b)## [1] 2
min(b)## [1] 1
sd(b)## [1] 0.7071068
var(b)## [1] 0.5
mean(b)## [1] 1.5
median(b)## [1] 1.5
quantile(b,.9)## 90%
## 1.9
sample(b,1) ## [1] 1
Erzeugen Sie einen Vektor b mit den Zahlen von 1 bis 5 und berechnen Sie…
den Mittelwert
die Varianz
die Standardabweichung
die quadratische Wurzel aus dem Mittelwert
| Datentyp | Beschreibung | Beispiel |
|---|---|---|
| numeric | ganze und reele Zahlen | 5, 3.462 |
| logical | logische Werte | FALSE, TRUE |
| character | Buchstaben und Zeichenfolgen | “Hallo” |
Quelle: R. Münnich und M. Knobelspieß (2007): Einführung in das statistische Programmpaket R
b <- c(1,2) # numeric
log <- c(T,F) # logical
char <-c("A","b") # character
fac <- as.factor(c(1,2)) # factorMit str() bekommt man den Objekttyp.
A1 <- c(1,2,3,4)
A1## [1] 1 2 3 4
A1[1]## [1] 1
A1[4]## [1] 4
A1[1:3]## [1] 1 2 3
A1[-4]## [1] 1 2 3
Beispieldaten generieren:
AGE <- c(20,35,48,12)
SEX <- c("m","w","w","m")Diese beiden Vektoren zu einem data.frame verbinden:
Daten <- data.frame(Alter=AGE,Geschlecht=SEX)Anzahl der Zeilen/Spalten herausfinden
nrow(Daten) # Zeilen## [1] 4
ncol(Daten) # Spalten## [1] 2
Indizieren eines dataframe:
AA <- 4:1
A2 <- cbind(A1,AA)
A2[1,1]## A1
## 1
A2[2,]## A1 AA
## 2 3
A2[,1]## [1] 1 2 3 4
A2[,1:2]## A1 AA
## [1,] 1 4
## [2,] 2 3
## [3,] 3 2
## [4,] 4 1
A <- matrix(seq(1,100), nrow = 4)
dim(A)## [1] 4 25
A3 <- array(1:8,c(2,2,2))
A3## , , 1
##
## [,1] [,2]
## [1,] 1 3
## [2,] 2 4
##
## , , 2
##
## [,1] [,2]
## [1,] 5 7
## [2,] 6 8
A3[,,2]## [,1] [,2]
## [1,] 5 7
## [2,] 6 8
A4 <- list(A1,1)
A4## [[1]]
## [1] 1 2 3 4
##
## [[2]]
## [1] 1
A4[[2]]## [1] 1
# Ist 1 größer als 2?
1>2## [1] FALSE
1<2## [1] TRUE
1==2## [1] FALSE
Diese Operatoren eignen sich gut um Datensätze einzuschränken
Daten## Alter Geschlecht
## 1 20 m
## 2 35 w
## 3 48 w
## 4 12 m
Daten[AGE>20,]## Alter Geschlecht
## 2 35 w
## 3 48 w
Daten[SEX=="w",]## Alter Geschlecht
## 2 35 w
## 3 48 w
# gleiches Ergebnis:
Daten[SEX!="m",]## Alter Geschlecht
## 2 35 w
## 3 48 w
# Ergebnis in ein Objekt speichern
subDat <- Daten[AGE>20,]# mehrere Bedingeungen können mit
# & verknüpft werden:
Daten[AGE>18 & SEX=="w",]## Alter Geschlecht
## 2 35 w
## 3 48 w
# Sequenz von 1 bis 10
1:10## [1] 1 2 3 4 5 6 7 8 9 10
Daten[1:3,]## Alter Geschlecht
## 1 20 m
## 2 35 w
## 3 48 w
seq(-2,8,by=1.5)## [1] -2.0 -0.5 1.0 2.5 4.0 5.5 7.0
a <-seq(3,12,length=12)
b <- seq(to=5,length=12,by=0.2)
d <- 1:10
d <- seq(1,10,1)
d <- seq(length=10,from=1,by=1)# wiederhole 1 10 mal
rep(1,10)## [1] 1 1 1 1 1 1 1 1 1 1
rep("A",10)## [1] "A" "A" "A" "A" "A" "A" "A" "A" "A" "A"
?pastepaste(1:4)## [1] "1" "2" "3" "4"
paste("A", 1:6, sep = "")## [1] "A1" "A2" "A3" "A4" "A5" "A6"
help.start()help(name)?meanexample(lm)R-project + Was ich schon immer wissen wollte
install.packages("lme4")
library(lme4)Die wichtigsten Pakete zur Visualisierung mit R:
Weitere interessante Pakete
Paket für den Import/Export - foreign
Paket - Latex und R (xtable) (xtable Galerie)
install.packages("devtools")
library(devtools)
install_github("hadley/ggplot2")Gehen Sie auf <cran.r-project.org> und suchen Sie in dem Bereich, wo die Pakete vorgestellt werden, nach Paketen,…
R unterstützt von Haus aus schon einige wichtige Formate:
read.csv()read.fwf()read.delim()So findet man heraus, in welchem Verzeichnis man sich gerade befindet
getwd()So kann man das Arbeitsverzeichnis ändern:
Man erzeugt ein Objekt in dem man den Pfad abspeichert:
main.path <- "C:/" # Beispiel für Windows
main.path <- "/users/Name/" # Beispiel für Mac
main.path <- "/home/user/" # Beispiel für LinuxUnd ändert dann den Pfad mit setwd()
setwd(main.path)Bei Windows ist es wichtig Slashs anstelle von Backslashs zu verwenden.
library(foreign) ist für den Import von fremden Datenformaten nötigread.csv() genutzt werden um die Daten einzulesen.read.csv2() wegen der Komma-Separierung braucht.library(foreign)
?read.csv
?read.csv2Zunächst muss das Arbeitsverzeichnis gesetzt werden, in dem sich die Daten befinden:
Dat <- read.csv("schuldaten_export.csv")Wenn es sich um Deutsche Daten handelt:
Dat <- read.csv2("schuldaten_export.csv")Dateien können auch direkt aus dem Internet geladen werden:
link<- "http://www.statistik.at/web_de/static/
mz_2013_sds_-_datensatz_080469.sav"
?read.spss
Dat <- read.spss(link,to.data.frame=T)rioinstall.packages("rio")library("rio")
x <- import("mtcars.csv")
y <- import("mtcars.rds")
z <- import("mtcars.dta")install.packages("Rz")
library(Rz)Laden Sie den Datensatz mit einer geeigneten Funktion in Ihre Console.
Finden Sie heraus, wieviele Beobachtungen und Variablen der Datensatz umfasst.
read.X() Funktionen stehen viele write.X() Funktionen zur Verfügung.RData)A <- c(1,2,3,4)
B <- c("A","B","C","D")
mydata <- data.frame(A,B)save(mydata, file="mydata.RData")write.csv(mydata,file="mydata.csv") library(xlsx)
write.xlsx(mydata,file="mydata.xlsx") library(foreign)
write.dta(mydata,file="mydata.dta") rio Paketlibrary("rio")
export(mtcars, "mtcars.csv")
export(mtcars, "mtcars.rds")
export(mtcars, "mtcars.dta")Quick R für das Exportieren von Daten:
Hilfe zum Export auf dem cran Server
Public-Use-File (PUF) Datei zur öffentlichen Nutzung - meist stark anonymisierte Daten (Beispiele: FDZ, Statistik Portal, Meine Region )
Scientific-Use-File (SUF) - Datei zur wissenschaftlichen Nutzung - anonymisierte Daten, die zu wissenschaftlichen Zwecken und zur Sekundäranalyse genutzt werden können.
On-Site-Nutzung - Arbeitsplätze für Gastwissenschaftler - Kontrollierte Datenfernverarbeitung
Auf dem Portal datahub.io sind sehr viele Beispieldatensätze in verschiedenen Formaten abrufbar.
Weitere Portale: OpenGov, okfn, enigma, Amazon Web Services (AWS)
Umweltdaten (National climatic data center)
library("FAOSTAT")Public Use File für Soziales in den USA Social security administration
National health and nutrition examination survey
library(survey)
data(nhanes)library(datasets)Beispiel Erdbeben Datensatz:
head(quakes)library(UScensus2010)WDI - World Development Indicators (World Bank) - Einführung in das Paket
library(WDI)WDIsearch('gdp')[1:10,]dat <- WDI(indicator='NY.GDP.PCAP.KD', country=c('MX','CA','US'), start=1960, end=2012)
head(dat)OpenStreetMap (OSM) ist ein kollaboratives Projekt um eine editierbare Weltkarte zu erzeugen.
library(osmar)
api <- osmsource_api()
library(ggmap)cityC <- geocode("Berlin",source="google")
bb <- center_bbox(cityC$lon,cityC$lat,1000, 1000)
uaBerlin <- get_osm(bb, source = api)Ausschnitte von OpenStreetMap für einzelne Städte (metro extracts)
Liste möglicher Datenquellen für räumliche Analysen (weltweit, Deutschland )
SALB - Administrative Grenzen
Kartendaten (openaprs)
library(twitteR)
library(streamR)library(mapdata)
data(worldHiresMapEnv)
map('worldHires', col=1:10)library(GDELTtools)
test.filter <- list(ActionGeo_ADM1Code=c("NI", "US"), ActionGeo_CountryCode="US")
test.results <- GetGDELT(start.date="1979-01-01", end.date="1979-12-31",
filter=test.filter)Mehr Daten hier
link1 <- "http://openflights.svn.sourceforge.net/viewvc/openflights/
openflights/data/airports.dat"
airport <- read.csv(link1, header = F)
link2 <- "http://openflights.svn.sourceforge.net/viewvc/openflights/
openflights/data/routes.dat"
route <- read.csv(link2, header = F)Hafen Daten (Natural earth data)
Data on airports and an example on the usage in R
link <- "http://www.fa-technik.adfc.de/code/opengeodb/DE9.tab"
info <- read.csv(link,sep="\t",header=F)Im base Paket sind die wichtigsten Streuungsmaße enthalten:
var()sd()min() und max()range()ab <- rnorm(100); var(ab)## [1] 1.072547
sd(ab); range(ab)## [1] 1.035638
## [1] -2.520023 2.572660
min(ab)## [1] -2.520023
max(ab)## [1] 2.57266
NAs vorhanden muss dies der Funktion mitgeteilt werdenab[10] <- NA
var(ab)## [1] NA
Bei fehlenden Werten muss ein weiteres Argument mitgegeben werden:
var(ab,na.rm=T)## [1] 1.078947
table()table() sind auch Kreuztabellierungen möglich indem zwei Variablen durch Komma getrennt werden: table(x,y) liefert Häufigkeiten von y für gegebene Ausprägungen von xx <- sample(1:10,100,replace=T)
table(x)## x
## 1 2 3 4 5 6 7 8 9 10
## 11 7 10 11 8 14 7 9 16 7
musician <- sample(c("yes","no"),100,replace=T)?tabletable(x)## x
## 1 2 3 4 5 6 7 8 9 10
## 11 7 10 11 8 14 7 9 16 7
table(x,musician)## musician
## x no yes
## 1 6 5
## 2 4 3
## 3 7 3
## 4 6 5
## 5 3 5
## 6 6 8
## 7 4 3
## 8 5 4
## 9 6 10
## 10 5 2
data(esoph)
table(esoph$agegp)##
## 25-34 35-44 45-54 55-64 65-74 75+
## 15 15 16 16 15 11
prop.table() liefert die relativen Häufigkeitentable() Funktion geschrieben erhält man die relativen Häufigkeiten bezogen auf alle ZellenDie Funktion prop.table()
table(esoph$agegp,esoph$alcgp)##
## 0-39g/day 40-79 80-119 120+
## 25-34 4 4 3 4
## 35-44 4 4 4 3
## 45-54 4 4 4 4
## 55-64 4 4 4 4
## 65-74 4 3 4 4
## 75+ 3 4 2 2
prop.table?prop.tableprop.table(table(esoph$agegp,
esoph$alcgp),1)##
## 0-39g/day 40-79 80-119 120+
## 25-34 0.2666667 0.2666667 0.2000000 0.2666667
## 35-44 0.2666667 0.2666667 0.2666667 0.2000000
## 45-54 0.2500000 0.2500000 0.2500000 0.2500000
## 55-64 0.2500000 0.2500000 0.2500000 0.2500000
## 65-74 0.2666667 0.2000000 0.2666667 0.2666667
## 75+ 0.2727273 0.3636364 0.1818182 0.1818182
aggregate() Funktion können Kennwerte für Untergruppen erstellt werdenaggregate(x,by,FUN) müssen mindestens drei Argumente übergeben werden:aggregate(state.x77,by=list(state.region),mean)## Group.1 Population Income Illiteracy Life Exp Murder HS Grad
## 1 Northeast 5495.111 4570.222 1.000000 71.26444 4.722222 53.96667
## 2 South 4208.125 4011.938 1.737500 69.70625 10.581250 44.34375
## 3 North Central 4803.000 4611.083 0.700000 71.76667 5.275000 54.51667
## 4 West 2915.308 4702.615 1.023077 71.23462 7.215385 62.00000
## Frost Area
## 1 132.7778 18141.00
## 2 64.6250 54605.12
## 3 138.8333 62652.00
## 4 102.1538 134463.00
x: ein oder mehrere Beobachtungsvektor(en) für den der Kennwert berechnet werden soll
by: eine oder mehrere bedingende Variable(n)
FUN: die Funktion welche den Kennwert berechnet (z.B. mean oder sd)
xtabs() in eine schöne zweidimensionale Tabelle überführt werdenApplyDat <- cbind(1:4,runif(4),rnorm(4))apply(ApplyDat,1,mean)## [1] 0.2174516 0.3037641 1.0847166 1.6806857
apply(ApplyDat,2,mean)## [1] 2.5000000 0.4733878 -0.5084242
apply(ApplyDat,1,var)## [1] 0.6143683 2.4650734 3.4903401 4.0587141
apply(ApplyDat,1,sd)## [1] 0.7838165 1.5700552 1.8682452 2.0146251
apply(ApplyDat,1,range)## [,1] [,2] [,3] [,4]
## [1,] -0.5676269 -1.098575 -0.7326386 0.3651438
## [2,] 1.0000000 2.000000 3.0000000 4.0000000
apply(ApplyDat,1,length)## [1] 3 3 3 3
Für margin=1 die Funktion mean auf die Reihen angewendet,
Für margin=2 die Funktion mean auf die Spalten angewendet,
Anstatt mean können auch andere Funktionen wie var, sd oder length verwendet werden.
ApplyDat <- data.frame(Income=rnorm(5,1400,200),
Sex=sample(c(1,2),5,replace=T))ApplyDat## Income Sex
## 1 1336.306 2
## 2 1231.769 2
## 3 1374.520 1
## 4 1430.820 1
## 5 1280.671 1
tapply(ApplyDat$Income,ApplyDat$Sex,mean)## 1 2
## 1362.004 1284.037
tapply(ApplyDat$Income,
ApplyDat$Sex,function(x)x)## $`1`
## [1] 1374.520 1430.820 1280.671
##
## $`2`
## [1] 1336.306 1231.769
Die Benutzung von apply, tapply, etc. (Artikel bei R-bloggers)
Erstellen Sie eine Matrix A mit 4 Zeilen und 25 Spalten, die die Werte 1 bis 100 enthält. Analog dazu erstellen Sie eine Matrix B mit 25 Zeilen und 4 Spalten, die die Werte 1 bis 100 enthält.
Berechnen Sie mittels dem apply()-Befehl den Mittelwert und die Varianz für jede Zeile von A bzw. B.
Berechnen Sie mittels dem apply()-Befehl den Mittelwert und die Varianz für jede Spalte von A bzw. B.
Standardisieren ist eine häuge Transformation von Daten; dafür wird der Mittelwert von der entsprechenden Zeile o der Spalte abgezogen und durch die entsprechende Standardab- weichung geteilt. Somit b esitzen die Daten einen Mittelwert von 0 und eine Standardab- weichung von 1. Standardisieren Sie die Spalten der Matrix A .
install.packages("ctv")
library("ctv")
install.views("Bayesian")library(mlmRev)
data(Chem97)Wir erstellen ein Histogramm der Variable gcsescore:
?histhist(Chem97$gcsescore)png, pdf oder jpegpng("Histogramm.png")
hist(Chem97$gcsescore)
dev.off()hist() plottet ein Histogramm der Datenhist() hat noch sehr viel mehr Argumente, die alle (sinnvolle) default values haben| Argument | Bedeutung | Beispiel |
|---|---|---|
| main | Überschrift | main=“Hallo Welt” |
| xlab | x-Achsenbeschriftung | xlab=“x-Werte” |
| ylab | y-Achsenbeschriftung | ylab=“y-Werte” |
| col | Farbe | col=“blue” |
hist(Chem97$gcsescore,col="blue",
main="Hallo Welt",ylab="y-Werte", xlab="x-Werte")Weitere Argumente:
?plot
# oder
?parbarplot() erzeugt aus einer Häufigkeitstabelle einen BarplottabScore <- table(Chem97$score)barplot(tabScore)barplot(tabScore)barplot(tabScore,col=rgb(0,0,1))barplot(tabScore,col=rgb(0,1,0))barplot(tabScore,col=rgb(1,0,0))barplot(tabScore,col=rgb(1,0,0,.3))boxplot()boxplot() muss mindestens ein Beobachtungsvektor übergeben werden?boxplotboxplot() ein sog. Formel-Objekt übergeben werdenboxplot(Chem97$gcsescore~Chem97$gender)Violinplot
# Beispieldaten erzeugen
x <- rnorm(100)
y <- rnorm(100)vioplotlibrary(vioplot)
plot(x, y, xlim=c(-5,5), ylim=c(-5,5))
vioplot(x, col="tomato", horizontal=TRUE, at=-4,
add=TRUE,lty=2, rectCol="gray")
vioplot(y, col="cyan", horizontal=FALSE, at=-4,
add=TRUE,lty=2)vioplot - Das Ergebnislibrary(beanplot)
par(mfrow = c(1,2))
boxplot(count~spray,data=InsectSprays,col="blue")
beanplot(count~spray,data=InsectSprays,col="orange")data(iris)
head(iris)## Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 1 5.1 3.5 1.4 0.2 setosa
## 2 4.9 3.0 1.4 0.2 setosa
## 3 4.7 3.2 1.3 0.2 setosa
## 4 4.6 3.1 1.5 0.2 setosa
## 5 5.0 3.6 1.4 0.2 setosa
## 6 5.4 3.9 1.7 0.4 setosa
petal length and width - Blütenblatt Länge und Breite
sepal length and width - Kelchblatt Länge und Breite
# Pearson Korrelationskoeffizient
cor(iris$Sepal.Length,iris$Petal.Length)## [1] 0.8717538
cor().pairs(iris[,1:4])library("psych")
pairs.panels(iris[1:4],bg=c("red","yellow","blue")
[iris$Species],pch=21,main="")# Pearson Korrelationskoeffizient
cor(iris[,1:4]) ## Sepal.Length Sepal.Width Petal.Length Petal.Width
## Sepal.Length 1.0000000 -0.1175698 0.8717538 0.8179411
## Sepal.Width -0.1175698 1.0000000 -0.4284401 -0.3661259
## Petal.Length 0.8717538 -0.4284401 1.0000000 0.9628654
## Petal.Width 0.8179411 -0.3661259 0.9628654 1.0000000
# Kendall's tau (Rangkorrelation)
cor(iris[,1:4], method = "kendall") ## Sepal.Length Sepal.Width Petal.Length Petal.Width
## Sepal.Length 1.00000000 -0.07699679 0.7185159 0.6553086
## Sepal.Width -0.07699679 1.00000000 -0.1859944 -0.1571257
## Petal.Length 0.71851593 -0.18599442 1.0000000 0.8068907
## Petal.Width 0.65530856 -0.15712566 0.8068907 1.0000000
# Spearman's rho (Rangkorrelation)
cor(iris[,1:4], method = "spearman") ## Sepal.Length Sepal.Width Petal.Length Petal.Width
## Sepal.Length 1.0000000 -0.1667777 0.8818981 0.8342888
## Sepal.Width -0.1667777 1.0000000 -0.3096351 -0.2890317
## Petal.Length 0.8818981 -0.3096351 1.0000000 0.9376668
## Petal.Width 0.8342888 -0.2890317 0.9376668 1.0000000
library("lattice")
library("AER")
data(BankWages)
levelplot(table(BankWages$education,BankWages$job))mosaicplot(~ Sex + Age + Survived,
data = Titanic, color = TRUE)Flächen werden entsprechend der Residuen eingefärbt:
mosaicplot(~ Sex + Age + Survived,
data = Titanic, shade = TRUE)Sachs - Angewandte Statistik mit R
It is designed to meet most typical graphics needs with minimal tuning, but can also be easily extended to handle most nonstandard requirements.
http://stat.ethz.ch/R-manual/R-devel/library/lattice/html/Lattice.html
library("lattice");library("mlmRev")
data(Chem97)
histogram(~ gcsescore, data = Chem97) histogram(~ gcsescore | factor(score),data = Chem97)densityplot(~ gcsescore | factor(score), Chem97,
groups=gender,plot.points=FALSE,auto.key=TRUE)bwplot(factor(score) ~ gcsescore | gender, Chem97)bwplot(gcsescore ~ gender | factor(score), Chem97,
layout = c(6, 1))barchart(yield ~ variety | site, data = barley,
groups = year, layout = c(1,6), stack = TRUE,
auto.key = list(space = "right"),
ylab = "Barley Yield (bushels/acre)",
scales = list(x = list(rot = 45)))densityplot( ~ height | voice.part, data = singer, layout = c(2, 4),
xlab = "Height (inches)", bw = 5)qq(gender ~ gcsescore | factor(score), Chem97,
f.value = ppoints(100), type = c("p", "g"), aspect = 1)xyplot(Sepal.Length + Sepal.Width ~ Petal.Length + Petal.Width | Species,
data = iris, scales = "free", layout = c(2, 2),
auto.key = list(x = .6, y = .7, corner = c(0, 0)))splom(~iris[1:4], groups = Species, data = iris)super.sym <- trellis.par.get("superpose.symbol")
splom(~iris[1:4], groups = Species, data = iris,
panel = panel.superpose,
key = list(title = "Three Varieties of Iris",
columns = 3,
points = list(pch = super.sym$pch[1:3],
col = super.sym$col[1:3]),
text = list(c("Setosa", "Versicolor", "Virginica"))))parallelplot(~iris[1:4] | Species, iris)data <- read.csv("oecd.csv",header = TRUE)Überprufen Sie die Dimension der OECD-Daten.
Berechnen Sie die Mittelwerte und Varianzen der einzelnen Variablen mit einem geeigneten apply Befehl.
In welchem Land waren die meisten Jugendlichen mindestens zweimal betrunken? Wie hoch ist der maximale Prozentsatz?
In welchem Land ist die Sterblichkeit am geringsten? Wie hoch ist sie in diesem Land?
Erstellen Sie einen neuen Datensatz, der aufsteigend nach dem Einkommen geordnet ist. Speichern Sie diesen in einer neuen .csv Datei
Maindonald - DataAnalysis
John H. Maindonald and W. John Braun
DAAG - Data Analysis and Graphics Data and Functions
install.packages("DAAG")library("DAAG")
data(roller)help on roller data:
?rollerSchätzen eines Regressionsmodells:
roller.lm <- lm(depression ~ weight, data = roller)So bekommt man die Schätzwerte:
summary(roller.lm)##
## Call:
## lm(formula = depression ~ weight, data = roller)
##
## Residuals:
## Min 1Q Median 3Q Max
## -8.180 -5.580 -1.346 5.920 8.020
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.0871 4.7543 -0.439 0.67227
## weight 2.6667 0.7002 3.808 0.00518 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 6.735 on 8 degrees of freedom
## Multiple R-squared: 0.6445, Adjusted R-squared: 0.6001
## F-statistic: 14.5 on 1 and 8 DF, p-value: 0.005175
Falls das Modell ohne Intercept geschätzt werden soll:
lm(depression ~ -1 + weight, data = roller)##
## Call:
## lm(formula = depression ~ -1 + weight, data = roller)
##
## Coefficients:
## weight
## 2.392
summary(roller.lm)##
## Call:
## lm(formula = depression ~ weight, data = roller)
##
## Residuals:
## Min 1Q Median 3Q Max
## -8.180 -5.580 -1.346 5.920 8.020
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.0871 4.7543 -0.439 0.67227
## weight 2.6667 0.7002 3.808 0.00518 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 6.735 on 8 degrees of freedom
## Multiple R-squared: 0.6445, Adjusted R-squared: 0.6001
## F-statistic: 14.5 on 1 and 8 DF, p-value: 0.005175
roller.lm ist nun ein spezielles Regressions-Objektpredict(roller.lm) # Vorhersage## 1 2 3 4 5 6 7
## 2.979669 6.179765 6.713114 10.713233 12.046606 14.180002 14.980026
## 8 9 10
## 18.180121 24.046962 30.980502
resid(roller.lm) # Residuen## 1 2 3 4 5 6
## -0.9796695 -5.1797646 -1.7131138 -5.7132327 7.9533944 5.8199976
## 7 8 9 10
## 8.0199738 -8.1801213 5.9530377 -5.9805017
plot(roller.lm,1)plot(roller.lm,2)Regression - r-bloggers
Das Komplette Buch von Faraway- sehr intuitiv geschrieben.
Gute Einführung auf Quick-R
Beschrieben wird Wegstrecke, dreier Spielzeugautos die in unterschiedlichen Winkeln Rampe herunterfuhren.
Lesen Sie den Datensatz toycars in einen dataframe data ein und wandeln Sie die Variable car des Datensatzes in einen Faktor (as.factor) um.
Erstellen Sie drei Boxplots, die die zurückgelegte Strecke getrennt nach dem Faktor car darstellen.
Schätzen Sie für jedes der 3 Autos separat die Parameter des folgenden linearen Mo dells mit Hilfe der Funktion lm()
\[ distance_i= \beta_0 + \beta_1 \cdot angle_i + \epsilon_i\]
Beispiele mit R-code
Faraway - Extending the linear model with r
Faraway - Practical Regression and Anova using R
glmglm()glm() muss 1. ein Formel-Objekt mitgegeben werden und 2. die Klasse (binomial, gaussian, Gamma) samt link-Funktion (logit, probit, cauchit, log, cloglog)install.packages("HSAUR")library("HSAUR")
data("plasma", package = "HSAUR")plasma_glm_1 <- glm(ESR ~ fibrinogen, data = plasma,
family = binomial())install.packages("faraway")library("faraway")data(orings)| temp | damage |
|---|---|
| 53 | 5 |
| 57 | 1 |
| 58 | 1 |
probitmod <- glm(cbind(damage,6-damage) ~ temp,
family=binomial(link=probit), orings)modp <- glm(Species ~ .,family=poisson,gala)MASS:library("MASS")
house.plr<-polr(Sat~Infl,weights=Freq,data=housing)data <- c(1,2,2,3,1,2,3,3,1,2,3,3,1)
fdata <- factor(data)
fdata## [1] 1 2 2 3 1 2 3 3 1 2 3 3 1
## Levels: 1 2 3
labels direkt definierenrdata <- factor(data,labels=c("I","II","III"))
rdata## [1] I II II III I II III III I II III III I
## Levels: I II III
levels(fdata) <- c('I','II','III')
fdata## [1] I II II III I II III III I II III III I
## Levels: I II III
mons <- c("March","April","January",
"November","January","September",
"October","September","November",
"August","January","November",
"November","February","May",
"August","July","December",
"August","August","September",
"November","February","April")
mons <- factor(mons)
table(mons)## mons
## April August December February January July March
## 2 4 1 2 3 1 1
## May November October September
## 1 5 1 3
mons <- factor(mons,levels=c("January","February",
"March","April","May","June",
"July","August","September",
"October","November",
"December"),
ordered=TRUE)
mons[1] < mons[2]## [1] TRUE
fert <- c(10,20,20,50,10,20,10,50,20)
fert <- factor(fert,levels=c(10,20,50),ordered=TRUE)
fert## [1] 10 20 20 50 10 20 10 50 20
## Levels: 10 < 20 < 50
mean(fert)## [1] NA
mean(as.numeric(levels(fert)[fert]))## [1] 23.33333
lets <- sample(letters,size=100,replace=TRUE)
lets <- factor(lets)
table(lets[1:5])##
## a b c d e f g h i j k l m n o p q r s t u v w x y z
## 0 1 1 0 0 0 0 0 0 0 0 1 0 1 0 0 0 1 0 0 0 0 0 0 0 0
ggplot2ggplot2<www.r-bloggers.com/basic-introduction-to-ggplot2/>
install.packages("ggplot2")library(ggplot2)diamonds Datensatzhead(diamonds)| carat | cut | color | clarity | depth | table | price | x | y | z |
|---|---|---|---|---|---|---|---|---|---|
| 0.23 | Ideal | E | SI2 | 61.5 | 55 | 326 | 3.95 | 3.98 | 2.43 |
| 0.21 | Premium | E | SI1 | 59.8 | 61 | 326 | 3.89 | 3.84 | 2.31 |
| 0.23 | Good | E | VS1 | 56.9 | 65 | 327 | 4.05 | 4.07 | 2.31 |
| 0.29 | Premium | I | VS2 | 62.4 | 58 | 334 | 4.20 | 4.23 | 2.63 |
| 0.31 | Good | J | SI2 | 63.3 | 58 | 335 | 4.34 | 4.35 | 2.75 |
| 0.24 | Very Good | J | VVS2 | 62.8 | 57 | 336 | 3.94 | 3.96 | 2.48 |
qplotqplot wird für schnelle Graphiken verwendet (quick plots)ggplot kann man alles bis ins Detail kontrollieren# histogram
qplot(depth, data=diamonds)qplot(cut, depth, data=diamonds)qplot(factor(cyl), data=mtcars,geom="bar")qplot(data=diamonds,x=cut,y=depth,geom="boxplot")# scatterplot
qplot(carat, depth, data=diamonds)qplot(carat, depth, data=diamonds,color=cut)myGG<-qplot(data=diamonds,x=carat,y=depth,color=carat)
myGG + stat_smooth(method="lm")qplot(factor(cyl), data=mtcars, geom="bar") +
coord_flip()ggplot(diamonds, aes(clarity, fill=cut)) + geom_bar()Es wird das Paket RColorBrewer verwendet um die Farbpalette zu ändern
install.packages("RColorBrewer")library(RColorBrewer)
myColors <- brewer.pal(5,"Accent")
names(myColors) <- levels(diamonds$cut)
colScale <- scale_colour_manual(name = "cut",
values = myColors)p <- ggplot(diamonds,aes(carat, depth,colour = cut)) +
geom_point()
p + colScaleggsave("Graphik.jpg")Arten von räumlichen Daten:
Das R-paket ggmap wird im folgenden genutzt um verschiedene Kartentypen darzustellen.
Mit qmap kann man eine schnelle Karte erzeugen.
ggmap:devtools::install_github("dkahle/ggmap")
devtools::install_github("hadley/ggplot2")
install.packages("ggmap")librarylibrary(ggmap)Und schon kann die erste Karte erstellt werden:
qmap("Mannheim")BBT <- qmap("Berlin Brandenburger Tor")
BBTqmap("Germany")qmap("Germany", zoom = 6)?qmapVerschiedene Abschnitte in der Hilfe:
Ausschnitt aus der Hilfe Seite zum Befehl qmap:
qmap Example
Das Beispiel kann man direkt in die Konsole kopieren:
# qmap("baylor university")
qmap("baylor university", zoom = 14)
# und so weiterqmap("Mannheim", zoom = 12)qmap('Mannheim', zoom = 13)qmap('Mannheim', zoom = 20)qmap('Mannheim', zoom = 14, maptype="satellite")qmap('Mannheim', zoom = 20, maptype="hybrid")qmap("Mannheim", zoom = 14, maptype="hybrid")Aus Physischen Karten kann man Informationen über Berge, Flüsse und Seen ablesen.
Farben werden oft genutzt um Höhenunterschiede zu visualisieren
qmap('Schriesheim', zoom = 14,maptype="terrain")New York
Abstraktion wird genutzt um nur die essentiellen Informationen einer Karte zu zeigen.
Bsp. U-Bahn Karten - wichtig sind Richtungen und wenig Infos zur Orientierung
Im folgenden werden Karten vorgestellt, die sich gut als Hintergrundkarten eignen.
qmap('Mannheim', zoom = 14,maptype="watercolor",source="stamen")qmap('Mannheim', zoom = 14,
maptype="toner",source="stamen")qmap('Mannheim', zoom = 14,
maptype="toner-lite",source="stamen")qmap('Mannheim', zoom = 14,
maptype="toner-hybrid",source="stamen")qmap('Mannheim', zoom = 14,
maptype="terrain-lines",source="stamen")RstudioExport
<- ist der Zuweisungspfeil um ein Objekt zu erzeugenMA_map <- qmap('Mannheim',
zoom = 14,
maptype="toner",
source="stamen")Geocoding (…) uses a description of a location, most typically a postal address or place name, to find geographic coordinates from spatial reference data …
library(ggmap)
geocode("Mannheim",source="google")| lon | lat |
|---|---|
| 8.463243 | 49.48604 |
| cities | lon | lat |
|---|---|---|
| Hamburg | 9.993682 | 53.55108 |
| Koeln | 6.960279 | 50.93753 |
| Dresden | 13.737262 | 51.05041 |
| Muenchen | 11.581981 | 48.13513 |
Reverse geocoding is the process of back (reverse) coding of a point location (latitude, longitude) to a readable address or place name. This permits the identification of nearby street addresses, places, and/or areal subdivisions such as neighbourhoods, county, state, or country.
Quelle: Wikipedia
revgeocode(c(48,8))## [1] "Unnamed Road, Somalia"
mapdist("Q1, 4 Mannheim","B2, 1 Mannheim")## from to m km miles seconds minutes
## 1 Q1, 4 Mannheim B2, 1 Mannheim 749 0.749 0.4654286 212 3.533333
## hours
## 1 0.05888889
mapdist("Q1, 4 Mannheim","B2, 1 Mannheim",mode="walking")## from to m km miles seconds minutes hours
## 1 Q1, 4 Mannheim B2, 1 Mannheim 546 0.546 0.3392844 423 7.05 0.1175
mapdist("Q1, 4 Mannheim","B2, 1 Mannheim",mode="bicycling")## from to m km miles seconds minutes
## 1 Q1, 4 Mannheim B2, 1 Mannheim 555 0.555 0.344877 215 3.583333
## hours
## 1 0.05972222
POI1 <- geocode("B2, 1 Mannheim",source="google")
POI2 <- geocode("Hbf Mannheim",source="google")
POI3 <- geocode("Mannheim, Friedrichsplatz",source="google")
ListPOI <-rbind(POI1,POI2,POI3)
POI1;POI2;POI3## lon lat
## 1 8.462844 49.48569
## lon lat
## 1 8.469879 49.47972
## lon lat
## 1 8.475208 49.48326
MA_map +
geom_point(aes(x = lon, y = lat),
data = ListPOI)MA_map +
geom_point(aes(x = lon, y = lat),col="red",
data = ListPOI)ListPOI$color <- c("A","B","C")
MA_map +
geom_point(aes(x = lon, y = lat,col=color),
data = ListPOI)ListPOI$size <- c(10,20,30)
MA_map +
geom_point(aes(x = lon, y = lat,col=color,size=size),
data = ListPOI)from <- "Mannheim Hbf"
to <- "Mannheim B2 , 1"
route_df <- route(from, to, structure = "route")qmap("Mannheim Hbf", zoom = 14) +
geom_path(
aes(x = lon, y = lat), colour = "red", size = 1.5,
data = route_df, lineend = "round"
)Wie fügt man Punkte hinzu
Nutzung von geom_point
Question on stackoverflow
pic
Was klar sein sollte:
visregEin einfaches Modell
N <- 5
x1 <- rnorm(N)
y <- runif(N)mod1 <- lm(y~x1)
pre <- predict(mod1)plot(x1,y)
abline(mod1)
segments(x1, y, x1, pre, col="red")visreg-PaketEin Modell wird auf dem airquality Datensatz geschätzt
install.packages("visreg")library(visreg)
fit <- lm(Ozone ~ Solar.R + Wind + Temp, data = airquality)visreg(fit)visreg visualisiert.visreg(fit, "Wind", type = "contrast")visregvisreg(fit, "Wind", type = "contrast")visreg-Paketvisreg(fit, "Wind", type = "conditional")Mit visreg können die Effekte bei Faktoren visualisiert werden.
airquality$Heat <- cut(airquality$Temp, 3,
labels=c("Cool", "Mild", "Hot"))
fit.heat <- lm(Ozone ~ Solar.R + Wind + Heat,
data = airquality)par(mfrow=c(1,2))
visreg(fit.heat, "Heat", type = "contrast")
visreg(fit.heat, "Heat", type = "conditional")airquality$Heat <- cut(airquality$Temp, 3,
labels=c("Cool", "Mild", "Hot"))
fit <- lm(Ozone ~ Solar.R + Wind * Heat, data = airquality)layoutvisreg(fit, "Wind", by = "Heat",layout=c(3,1))visreg - Interaktionen overlayfit <- lm(Ozone ~ Solar.R + Wind * Heat, data = airquality)
visreg(fit, "Wind", by="Heat", overlay=TRUE, partial=FALSE)visreg - visreg2dfit2 <- lm(Ozone ~ Solar.R + Wind * Temp, data = airquality)
visreg2d(fit2, "Wind", "Temp", plot.type = "image")visreg2d(fit2, "Wind", "Temp", plot.type = "persp")library("Rcmdr")Das Paket iplots
install.packages("iplots",dep=TRUE)library(iplots)cyl.f <- factor(mtcars$cyl)
gear.f <- factor(mtcars$factor)
attach(mtcars)
ihist(mpg) # histogram
ibar(carb) # barchart
iplot(mpg, wt) # scatter plot
ibox(mtcars[c("qsec","disp","hp")]) # boxplots
ipcp(mtcars[c("mpg","wt","hp")]) # parallel coordinates
imosaic(cyl.f,gear.f) # mosaic plot library(rggobi)
g <- ggobi(mydata) attach(mydata)
plot(x, y) # scatterplot
identify(x, y, labels=row.names(mydata)) # identify points
coords <- locator(type="l") # add lines
coords # display list library(stargazer)
stargazer(attitude)library(knitr)
kable(head(iris), format = "latex")pic
Folien zum Workshop:
https://github.com/Japhilko/npRegression/tree/master/slides
library(splines)interaktiven Karte und Rcode um eine interaktive Karte mit leaflet zu erzeugen.